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Abstract 

We propose a prescription based on the Fokker-Planck equation in the 
Stratonovich approach, with the diffusion coefficient dependent on tempo- 
ral and spatial coordinates, for describing heat conduction by phonons in 
small structures. This equation can be analytically solved for a broad class of 
diffusion coefficients. It can also describe non-Gaussian processes. Further, 
it generalizes the model investigated by Naqvi and Waldenstri;^m (PRL, 95 
(2005), 065901). We show that our solutions can fit well the results derived 
from the Boltzmann equation. 



1 



The processes of the heat conduction in physical systems involve the microscopic trans- 
ports of the energy carriers. In general, these processes can be demarcated by certain 
characteristic time and length scales of the energy carriers such as collision time, mean free 
time, relaxation time, diffusion time, mean free path, relaxation length and diffusion length. 
However, the descriptions of the heat conduction processes have encompassed several theo- 
retical approaches. The simplest one is based on the Fourier law [1]. This law enjoys some 
universality in the description of the heat transfer due to the fact that it can be applied 
to a wide range of the physical systems with different energy carriers, in the macroscales. 
The Fourier law breaks down for anomalous heat conduction systems, where the thermal 
conductivity k diverges with the system size L as L^. Theoretical investigation based on a 
connection between the anomalous diffusion processes and anomalous heat conduction in a 
one-dimensional systems has been carried out [2]. Also, the Fourier law breaks down in the 
domain of the microscales such as the heat transport in a thin film [3]. For this last case, 
the equation of phonon radiative transfer (EPRT) can describe well, however, it is difficult 
to be solved. Recently, several authors have attempted to replace easier models from which 
can give a good approximation of the EPRT results [4-6]. The ballistic-diffusive equations 
(BDE), derived from the Boltzmann equation under the relaxation time approximation, 
which has been introduced by Chen [4,5] can capture the behaviors of the temperature and 
heat flux of the EPRT. An other description is based on the Brownian motion which has been 
introduced by Naqvi and Waldenstrom (NW) [6] . This last model can describe the temper- 
ature of the EPRT very well, but its heat flux deviates visibly in some spatial range (see 
Fig. l),where t* = t/r, ^ = x/L, AT = Ti - To, 6 = {T - To) /AT , = q/lCvAT), r = l/v 
is the mean-free time, v is the average of sound, / is the mean-free path and Kn = l/L is 
the Knudsen number. 

In this letter, we propose to generahze the NW model. Our model is based on the Fokker- 
Planck equation in the Stratonovich approach [7] with the diffusion coefficient that depends 
on time and space. This model can be analytically solved for a broad class of diffusion 
coefficients, and it also presents interesting asymptotic properties [8]. We show that our 
solutions can give a good approximation of the EPRT results. 

In order to motivate our proposal, we first discuss the NW model which is given by 

dtT{x,t) = a{t)d^^T{x,t)-b{t)d^T{x,t) , (1) 

where T{x, t) is the temperature. In their analysis, Naqvi and Waldenstrom have chosen 
a{t*) = k{1 — e^**) and b{t) — 0, where k = vl/3 is the thermal diffusivity. In this case, 
the heat ffux is given by g = —{\/K)a{f*)dxT{x.t), where A = Cvl/3 and C is the specific 
heat per unit volume. For t* ^ 1, the model recovers the Fourier equation. As have 
been noted by the authors, the model (1) describes a Gaussian process and it can only 
approximate to the results of EPRT. For b{t) — 0, we can show that the solutions of Eq. 
(1) can not be improved for any choice of a(t), and it can only obtain the similar results 
described by a{t*) = k(1 — e^**). In order to show this fact, we plot A = —(j)/{d^6) against 
the nondimensional coordinate ^ (Fig.2) from the data of EPRT (Fig.l). We see that A 
does not remain constant. From the NW model we obtain A = Kn{l — e~^)/3 = 0.21. 
This value approximates to the first part of the curve of Fig.2 well. Then, around the value 
^ = 0.6, the curve begins to deviate visibly from the value 0.21. As can also be seen from 
Fig. 1 the heat fiux begins to deviate from the EPRT result around the value ^ = 0.6. This 
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shows that the model (1), for b{t) = 0, can not be improved. On the other hand, for any 
function of b{t) different from a{t), the solution of equation (1) may not be easily obtained 
and the method of separation of variables can not be used, either. 
In our proposal we consider the following equation 

dtT{x, t) = Ka{t*)d, [D{x)d, {D{x)T{x, t))] (2) 

and the heat flux given by 

q = -Xa{t*)D{x)d:c {D{x)T{x,t)) . (3) 

We note that if T{x, t) is replaced by the probability density, then Eq. (2) becomes a 
stochastic equation namely the Fokker-Planck equation in the Stratonovich approach which 
is obtained from the Langevin equation with a multiplicative noise term [7] . We see that Eq. 
(2) generahzes the NW model with b{t) = 0. In this case, we recover Eq. (1) for constant 
D{x). 

The application of Eqs. (2) and (3) to the physical systems is to consider a slab of thick- 
ness L coupled to two thermal reservoirs. At time t = 0, one face (at x = L) is maintained at 
the temperature Tq, whereas the other face (at a; = 0) is raised to a temperature Ti. More- 
over, initially the slab is maintained at a uniform temperature Tq. In order to compare with 
the results of other approaches we use the nondimensional quantities defined above. We also 
consider the diffusion coefficient only depends on the nondimensional variables. The solution 
of Eq. (2) in terms of the nondimensional variables can be obtained by the transformations: 
du/d^ = l/D{^) and ds/dt* = a{t*). Then, Eq. (2) reduces to 9,p(C, s) = {KT/L^)dlp{u,s), 
where p — DO. For convenience, we set m(^ = 0) = and s(r = 0) = 0; and the solutions 
for 9{^,t*) and 0(^,r) subject to the boundary conditions above are given by 



and 



(4) 



where Dq is the value of D{^) at ^ = and Ui is the value of u{^) at ^ = 1. For D = 1, we 
recover the results of NW model 

r 2 - sm(/»^Oexpf- "'"'^"':^'^^'^^ )l 

^(e, = 1 - c - ^ E 1 — ^ (6) 



and 

t*) = ^ |l + 2 E cos (mTrO exp ^ 3' ' ' j | ■ (7) 

For our numerical investigation we choose a{t*) — k{1 — e~^^*) and D{^) = 
Pi (1 / (1 where h, pi, p2 and ps are the parameters to be adjusted. For 
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simplicity, we have investigated the numerical solutions by using n as an integer. It seems 
that the numerical results, for n = 4, are better than other values of n. In Fig. 3 and Fig. 
4, we show the temperature and heat flux distributions obtained from Eqs. (4) and (5) by 
using Di{C) = (1 +P2C^)/(1 +P3C^) and £'2(0 = Pa (1 +P2C^) /b2 (1 +^3^^)], respectively. 
Both the results can fit the EPRT results well. We see that our numerical solutions, by using 
Di, can fit the results of EPRT better than those obtained by D2. However, the advantage 
of D2 is that it tends to unity for ,^ ^ 1. As have been noted by Joshi and Majumdar [3] , 
when the size of a slab is much larger than the phonon mean free path the heat transport 
can be modeled by the Fourier law. 

In Fig. 5 we compare the temporal behaviors of the nondimensional heat flux (at ^ = 0) 
obtained from the BME and FPE. We note that the EPRT data have not been included in 
Fig. 5 due to the fact that they have not been rescaled [6]. The upper curves (Kn = 10) 
show that our result has a similar behavior of that described by the BME, but our curve is 
lower than that described by the BME for not too small t*. The middle curves {Kn — 1) 
also show that our curve is lower than that of the BME for t* > 1. We note that the BME 
curve has a pronounced minimum around the value t* — 1.8, whereas our curve decays to a 
lower point and then keep it nearly straight. As a consequence, our curve can qualitatively 
reproduce the behavior of the EPRT result (see Fig. 2 of Ref. [4]) better than that of the 
BME. Finally, the lower curves show that our curve reproduces the BME result. 

As have been demonstrated in [3,4,6] the Fourier equation and Cattaneo equation [9] can 
not describe the EPRT results and consequently they fail to describe the heat conduction 
on small scales such as in a thin film. Moreover, the Fourier equation leads to a divergent 
heat flux for t* ^ and the Cattaneo equation produces artificial heat flux oscillation. 
In this work we have mainly concentrated our studies on the MBE and FPE, and their 
results have been compared with those calculated by using the EPRT and BDE. Our results 
based on the FPE can describe well the EPRT results. We have chosen the FPE in the 
Stratonovich approach due to the fact that it can be analytically solved for a broad class of 
diffusion coefficients, whereas in other approaches such as Ito and postpoint discretization 
approaches (see [8] and the references therein) we do not have the same facility. However, 
there is no reason to choose solely the FPE in the Stratonovich approach. In fact, numerical 
calculation can be used to obtain the solutions of other approaches. In order to choose 
which of the FPE approaches is more adequate for describing the heat conduction on small 
scales we need further information of the microscopic structure of the systems. We note that 
the FPE with the diffusion coefficient that depends on time and space can describe non- 
Gaussian processes with white noise [8]. In particular, the non-Gaussian processes obtained 
from the FPE are due to the dependence of the diffusion coefficient on the spatial coordinate. 
Therefore, Eq. (1) can only describe Gaussian processes and it is a particular case of Eq. 
(2). Finally, we would mention that we have used the coefficients a{t*) = (1 — e"'** ) and 
-0(0 = Pi (1 +P2C") / (1 + PsC") because they are simple expressions and they can fit well 
the EPRT results. However, other more elaborate expressions may also be employed for 
improving our results above. 
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FIGURE CAPTIONS 



Fig. 1 - Behaviors of the nondimensional temperature and heat flux of Brownian motion 
equation (MBE), Equation of phonon radiative transfer (EPRT) and baUistic diffusive equa- 
tions (BDE) for = 1 and = 1. The data of the EPRT and BDE have been extracted 
from Ref. [4]. 

Fig. 2 - Plot of the ratio A = —0/(9^6') in function of the nondimensional coordinate ^ 
calculated from the data of EPRT (Fig.l). 

Fig. 3 - Comparison of the nondimensional temperature and heat flux in terms of the 
nondimensional coordinate ^ obtained from the EPRT, FPE and BDE for Kn = 1 and 
t* — 1. In the case of FPE, the results have been calculated by using a{t*) = (1 — e"'***) 
and — (1 +^2^^) / (1 +^3^^) with the parameters given hy h — 1.01, p2 — 0.7 and 

Pa = 2.4. 

Fig. 4 - Comparison of the nondimensional temperature and heat flux in terms of the 
nondimensional coordinate ^ obtained from the EPRT, FPE and BDE for Kn = 1 and 

t* = 1. In the case of FPE, the results have been calculated by using a{t*) = (1 — e~'***) and 
-^2(0 = P3 (1 +P2^^) /\p2 (1 +^3^^)] with the parameters given hy h — 0.265, p2 — 0.63 and 
P3 = 1-015. 

Fig. 5 - Comparison of the nondimensional heat flux in function of the nondimensional 
time f* obtained from the BME and FPE at = 0. In this plot we have used a(f*) = 
(1 — e~^**) and -Di(C) = (1 +^2^"^) / (1 +^3^"^) with the parameters given hy h = 1.01, 
P2 = 0.7 and ps — 2.4. The sohd hues correspond to the BME data, whereas the dotted lines 
correspond to the FPE data. The pairs of the curves from top to bottom are calculated by 
using Kn — 10, 1, 0.1, respectively. 
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